Synergistic Charge Transfer Effect in Ferrous Heme–CO Bonding within Cytochrome P450

We conducted ab initio valence bond (VB) calculations employing the valence bond self-consistent field (VBSCF) and breathing orbital valence bond (BOVB) methods to investigate the nature of the coordination bonding between ferrous heme and carbon monoxide (CO) within cytochrome P450. These calculations revealed the significant influence exerted by both proximal and equatorial ligands on the π-backdonation effect from the heme to the CO. Moreover, our VB calculations unveiled a phenomenon of synergistic charge transfer (sCT). In the case of ferrous heme–CO bonding, the significant stabilization in this sCT arises from cooperative resonance between the VB structures associated with σ donation and π backdonation. Unlike many other ligands, CO possesses the unique ability to establish two mutually perpendicular π-backdonation orbital interaction pairs, leading to an intensified stabilization attributed to σ–π resonance. Furthermore, while of a smaller energy magnitude, sCT due to one π–π pair is also present, contributing to the differential stabilization of ferrous heme–CO bonding.


Introduction
The interaction between the ferrous heme and carbon monoxide (CO) stands as one of the most fundamental types of coordination bonds found in cytochrome P450 enzymes (P450s).The distinctive spectral features displayed by the ferrous-CO complex of P450s not only served as the inspiration behind the unique nomenclature "P450" but also played a crucial role in the initial discovery and characterization of these enzymes [1][2][3].The nature of ferrous heme-CO bonding in P450s and other heme proteins has been extensively studied using techniques such as infrared and resonance Raman spectroscopy [4,5].Furthermore, the strong affinity of ferrous heme for CO remains a widely utilized phenomenon in contemporary experimental investigations, including CO binding assays and deliberate enzyme inhibition strategies [6][7][8].
At the core of ferrous heme-CO bonding lies the pivotal role played by bidirectional charge transfer (CT) [4,5,[9][10][11].Among the array of ligands capable of binding to the P450 heme, CO distinguishes itself due to its ability to establish two mutually perpendicular pairs of orbital interactions for π-backdonation CT, as depicted in Figure 1.To quantitatively assess the effects of σ donation and π backdonation, as well as to determine their relative significance-unattainable through spectroscopic studies-we have recently conducted comprehensive theoretical studies employing complementary occupied-virtual orbital pair (COVP) and quantum theory of atoms-in-molecules (QTAIM) analyses based on density functional theory (DFT) [10,12,13].Through these analyses, we have convincingly demonstrated that π backdonation holds greater significance than σ donation in the ferrous heme-CO complex within P450s.Furthermore, our more recent research has shown that the energetic stabilization for π backdonation in the ferrous heme-CO complex surpasses what could be expected from the amount of transferred electron density [11].This intriguing observation suggests the presence of additional stabilization in the ferrous heme-CO complex when compared to other heme-ligand complexes.Thus, the nature of ferrous heme-CO bonding in P450s presents an intellectually captivating subject, worthy of renewed scrutiny.
Molecules 2024, 29, x FOR PEER REVIEW 2 of 13 complex surpasses what could be expected from the amount of transferred electron density [11].This intriguing observation suggests the presence of additional stabilization in the ferrous heme-CO complex when compared to other heme-ligand complexes.Thus, the nature of ferrous heme-CO bonding in P450s presents an intellectually captivating subject, worthy of renewed scrutiny.In our quest to unravel the quantum mechanical foundations of ferrous heme-CO coordination bonding in P450s, this study leverages ab initio valence bond (VB) theory.Beyond its indisputable role in advancing our comprehension of chemistry after the establishment of quantum mechanics [14], VB theory has undergone significant development in recent decades, becoming a quantitatively competitive ab initio tool [15][16][17].As it has evolved, ab initio VB theory has extended its application beyond organic molecules, providing invaluable insights into donor-acceptor and metal-ligand interactions [18][19][20][21][22][23][24][25].Furthermore, VB calculations have brought to light a novel phenomenon termed chargeshift bonding, significantly enriching our comprehension of chemical bonds [26,27].The driving force behind charge-shift bonding defies conventional covalent bonding models, with resonance between covalent and ionic VB structures emerging as a prominent factor facilitating bonding.The potential and utility of ab initio VB theory are poised for further expansion into the realm of complex biological systems containing transition metals, such as P450s.
This study pursues two primary objectives: firstly, to evaluate how proximal and equatorial ligands influence the relative significance of σ donation and π backdonation within the framework of bidirectional CT in P450s; and secondly, to elucidate the underlying mechanisms responsible for the observed additional CT stabilization in the ferrous heme-CO complex.

Wave Function Analysis: Relative Significance of Different VB Structures
In our VB analysis of ferrous heme-CO bonding, we employed the three models (I-III) depicted in Scheme 1.Further details regarding the model setup can be found in the Materials and Methods section (Section 3). Figure 2 displays the hybrid atomic orbitals (AOs) associated with the σ-donation and π-backdonation bonding pairs in model III, as determined through the VB(all) calculation with the valence bond self-consistent field (VBSCF) method.Similar orbitals were obtained for models I and II (Figures S2 and S3).These orbitals align with the conventional understanding of σ donation occurring from the n(CO) orbital to the  (Fe) orbital and π backdonation from the  (Fe) and  (Fe) orbitals to the π*(CO) orbitals (Figure 1).In our quest to unravel the quantum mechanical foundations of ferrous heme-CO coordination bonding in P450s, this study leverages ab initio valence bond (VB) theory.Beyond its indisputable role in advancing our comprehension of chemistry after the establishment of quantum mechanics [14], VB theory has undergone significant development in recent decades, becoming a quantitatively competitive ab initio tool [15][16][17].As it has evolved, ab initio VB theory has extended its application beyond organic molecules, providing invaluable insights into donor-acceptor and metal-ligand interactions [18][19][20][21][22][23][24][25].Furthermore, VB calculations have brought to light a novel phenomenon termed charge-shift bonding, significantly enriching our comprehension of chemical bonds [26,27].The driving force behind charge-shift bonding defies conventional covalent bonding models, with resonance between covalent and ionic VB structures emerging as a prominent factor facilitating bonding.The potential and utility of ab initio VB theory are poised for further expansion into the realm of complex biological systems containing transition metals, such as P450s.
This study pursues two primary objectives: firstly, to evaluate how proximal and equatorial ligands influence the relative significance of σ donation and π backdonation within the framework of bidirectional CT in P450s; and secondly, to elucidate the underlying mechanisms responsible for the observed additional CT stabilization in the ferrous heme-CO complex.

Wave Function Analysis: Relative Significance of Different VB Structures
In our VB analysis of ferrous heme-CO bonding, we employed the three models (I-III) depicted in Scheme 1.Further details regarding the model setup can be found in the Materials and Methods section (Section 3). Figure 2 displays the hybrid atomic orbitals (AOs) associated with the σ-donation and π-backdonation bonding pairs in model III, as determined through the VB(all) calculation with the valence bond self-consistent field (VBSCF) method.Similar orbitals were obtained for models I and II (Figures S2 and S3).These orbitals align with the conventional understanding of σ donation occurring from the n(CO) orbital to the d z 2 (Fe) orbital and π backdonation from the d xz (Fe) and d yz (Fe) orbitals to the π*(CO) orbitals (Figure 1).Table 1 presents the weights (WK) of the representative VB structures illustrated in Scheme 2 (see also Section 3 for theoretical details).The weights were calculated from the VB structure coefficients and the overlap integrals between VB structures, using the Coulson-Chirgwin formula [28].A noteworthy observation from the table is that VBSCF-calculated W1 consistently demonstrates the largest weight of Φ1 (68.39% in model I, 69.26% in model II, and 50.25% in model III) (Table 1a).In model I, W2 is the second largest and significantly larger than the W3 and W4 values.These results suggest that σ donation prevails over π backdonation in this model.In model II, W2 decreases while W3 and W4 increase compared to their corresponding weights in model I, indicating that the proximal HS -ligand attenuates σ donation while enhancing π backdonation.The attenuation of σ Scheme 1. Models I-III for VB analysis, constructed utilizing the previously optimized geometry of a P450 model [10].Table 1 presents the weights (WK) of the representative VB structures illustrated in Scheme 2 (see also Section 3 for theoretical details).The weights were calculated from the VB structure coefficients and the overlap integrals between VB structures, using the Coulson-Chirgwin formula [28].A noteworthy observation from the table is that VBSCF-calculated W1 consistently demonstrates the largest weight of Φ1 (68.39% in model I, 69.26% in model II, and 50.25% in model III) (Table 1a).In model I, W2 is the second largest and significantly larger than the W3 and W4 values.These results suggest that σ donation prevails over π backdonation in this model.In model II, W2 decreases while W3 and W4 increase compared to their corresponding weights in model I, indicating that the proximal HS -ligand attenuates σ donation while enhancing π backdonation.The attenuation of σ Table 1 presents the weights (W K ) of the representative VB structures illustrated in Scheme 2 (see also Section 3 for theoretical details).The weights were calculated from the VB structure coefficients and the overlap integrals between VB structures, using the Coulson-Chirgwin formula [28].A noteworthy observation from the table is that VBSCFcalculated W 1 consistently demonstrates the largest weight of Φ 1 (68.39% in model I, 69.26% in model II, and 50.25% in model III) (Table 1a).In model I, W 2 is the second largest and significantly larger than the W 3 and W 4 values.These results suggest that σ donation prevails over π backdonation in this model.In model II, W 2 decreases while W 3 and W 4 increase compared to their corresponding weights in model I, indicating that the proximal HS -ligand attenuates σ donation while enhancing π backdonation.The attenuation of σ donation and the enhancement of π backdonation are even more pronounced in model III, where W 2 -W 4 exhibit similar values (8.59, 8.52, and 8.53%, respectively), underscoring the role of equatorial ligands in promoting π backdonation.Another intriguing finding for model III is the relatively large values of W 5 and W 6 (9.27%), indicating that states concurrently involving one Heitler-London (HL) bonding pair for σ donation and one for π backdonation emerge as significant contributors to the overall wave function.Breathing orbital valence bond (BOVB) calculations produced results qualitatively similar to those obtained from VBSCF calculations (Table 1b).However, the BOVB decreased the weights of Φ 1 while increasing the weights of Φ 2 -Φ 4 , suggesting that the orbital breathing effect enhances CT in both directions.W 7 and W 8 were found to be very small in both the VBSCF and BOVB calculations.1b).However, the BOVB decreased the weights of Φ1 while increasing the weights of Φ2-Φ4, suggesting that the orbital breathing effect enhances CT in both directions.W7 and W8 were found to be very small in both the VBSCF and BOVB calculations.To assess the relative significance of σ donation and π backdonation, we introduced metrics to evaluate their effective weights, denoted as W(σ) and W(π), respectively.These metrics are defined as follows: To assess the relative significance of σ donation and π backdonation, we introduced metrics to evaluate their effective weights, denoted as W(σ) and W(π), respectively.These metrics are defined as follows: Table 2 provides a summary of the calculated W(σ) and W(π) values for all models.These values reaffirm the dominance of σ donation in model I (W(σ)>>W(π)), while models II and III exhibit an increased relative significance of π backdonation.Notably, in the case of model III, W(π) surpasses W(σ), which is consistent with our previous finding that π backdonation dominates over σ donation in ferrous heme-CO bonding in P450s [10,11].
Thus, the proximal and equatorial ligands exert substantial influence on amplifying π backdonation.However, it should also be noted that despite the decrease in the W 2 value when transitioning from model II to model III (Table 1), the W(σ) value for model III surpasses that of model II (Table 2).This outcome stems from the heightened contribution of the VB structures simultaneously involving both σ and π bonds, specifically Φ 5 and Φ 6 .Once again, the BOVB calculations yielded qualitatively similar results to the VBSCF calculations, albeit with elevated values for W(σ) and W(π).

Energy Analysis: Resonance between VB Structures
To gain further insights into the resonance stabilization resulting from the mixing of different VB structures, we conducted energy analyses.Resonance, in essence, involves the mixing of distinct VB structures (diabatic states) that represent different bonding patterns, resulting in a stabilizing effect on the overall wave function.For example, the remarkable stability of benzene finds its rationale in the resonance between the two Kekulé VB structures (alongside other minor ones) [15].Specifically, we calculated the resonance energy (RE) in various states (X) with respect to Φ 1 , the most dominant VB structure, using the following equation: In this equation, the energy of the variationally optimized Φ 1 structure, E(Φ 1 ), is calculated as Φ 1 Ĥ |Φ 1 using the Hamiltonian Ĥ.The energy of other states, E(X), is also determined in the same manner.For example, X in this equation is σπ x π y in the case of VB(all) calculations, but other states from deactivated VB calculations, such as VB(σ), can also be used as X.It is important to note that a more positive value of RE indicates a larger resonance stabilization effect.Table 3a provides a summary of the VBSCF-calculated RE values for different states.When we compare the RE(σ), RE(π x ), and RE(π y ) values for model I, which gauge the energetic stabilization attributable to respective CT processes, it becomes evident that σ donation plays a notably more substantial role in stabilizing the system compared to π backdonation.However, in model II, the relative significance of σ donation and π backdonation undergoes a notable change.Here, the RE(σ) value decreases from 26.98 kcal/mol (model I) to 17.49 kcal/mol (model II), whereas the RE(π x ) and RE(π y ) values increase from 4.03 kcal/mol (model I) to 7.84-7.86kcal/mol (model II).This trend becomes even more pronounced in the case of model III, which exhibits RE(σ), RE(π x ), and RE(π y ) values of 15.77, 12.91, and 12.94 kcal/mol, respectively.Thus, the energetic stabilization for each of these CT routes in model III is of similar magnitude.These results highlight the critical roles played by proximal and equatorial ligands in energetically modulating both the σ-donation and π-backdonation effects, which can be attributed to alterations in the stability of the involved d-orbitals within distinct ligand fields [29].
In Figure 3, we present a comparison of the energy levels of ten molecular orbitals (MOs) in the frontier region (HOMO−4 to LUMO+4) obtained through B3LYP-D3BJ/def2-TZVP(6D,10F) calculations for models I-III after removing their CO ligand (additional details about these MOs can be found in Figure S5).As can be seen from this figure, d-orbitals are destabilized to a greater extent when more ligands coordinate to the ferrous center.This results in a higher electron-donating capacity of the iron for π backdonation and a reduced electron-accepting ability for σ donation.Real P450s feature a porphyrin dianion ligand rather than neutral NH 3 ligands, which would further accentuate the equatorial ligand effect on π backdonation.In fact, our earlier DFT study indicated that in a porphinebased P450 model, π backdonation was approximately twice as significant as σ donation in terms of energetic stabilization [10,11], instead of the approximately 1.6-fold ratio observed in model III.In Figure 3, we present a comparison of the energy levels of ten molecular orbitals (MOs) in the frontier region (HOMO−4 to LUMO+4) obtained through B3LYP-D3BJ/def2-TZVP(6D,10F) calculations for models I-III after removing their CO ligand (additional details about these MOs can be found in Figure S5).As can be seen from this figure, dorbitals are destabilized to a greater extent when more ligands coordinate to the ferrous center.This results in a higher electron-donating capacity of the iron for π backdonation and a reduced electron-accepting ability for σ donation.Real P450s feature a porphyrin dianion ligand rather than neutral NH3 ligands, which would further accentuate the equatorial ligand effect on π backdonation.In fact, our earlier DFT study indicated that in a porphine-based P450 model, π backdonation was approximately twice as significant as σ donation in terms of energetic stabilization [10,11], instead of the approximately 1.6-fold ratio observed in model III.

The data in Table 3 further reveal that the RE(σπx), RE(σπy), and RE(πxπy) values consistently surpass the corresponding RE values for one of their constituent CT routes, i.e., RE(σ), RE(πx), or RE(πy). In addition, the RE(σπxπy) value exceeds the values of RE(σπx), RE(σπy), and RE(πxπy).
These findings illustrate the stabilizing influence that results from the mixing of a larger number of VB structures.The data in Table 3 further reveal that the RE(σπ x ), RE(σπ y ), and RE(π x π y ) values consistently surpass the corresponding RE values for one of their constituent CT routes, i.e., RE(σ), RE(π x ), or RE(π y ).In addition, the RE(σπ x π y ) value exceeds the values of RE(σπ x ), RE(σπ y ), and RE(π x π y ).These findings illustrate the stabilizing influence that results from the mixing of a larger number of VB structures.
E(Φ 1 ) is identical in our VBSCF and BOVB calculations, and the variational space is broader in BOVB for other states.Hence, the RE values calculated using BOVB necessarily exceeded their corresponding VBSCF counterparts, as shown in Table 3b.In essence, this trend signifies that the VB structures involving CT experience greater stabilization in the BOVB framework as a result of the orbital breathing effect.Overall, while the enhancement of the RE values in BOVB compared to the VBSCF values is not particularly substantial, relatively significant increases are observed when π x and π y VB structures were simultaneously included in the wave function in model III.It is also noteworthy that both VBSCF and BOVB yield similar enhancements or reductions in RE values when transitioning from model I to model III.For instance, in the case of RE(π x ), the enhancement factor is approximately 3.2 in both VBSCF and BOVB results.
To obtain further insights into the resonance stabilization arising from the interaction between different VB structures associated with different CT modes, we conducted supplementary analyses by introducing ∆RE values for different combinations of VB structures (Equations ( 4)-( 6)).For instance, ∆RE(σπ x ), as defined in Equation ( 4), assesses the resonance stabilization arising from the combination of the σ and π x VB structures.
A more positive value of ∆RE again indicates a greater resonance-stabilizing effect.The calculated ∆RE values are presented in Figure 4, revealing a discernible pattern: the σ-π resonance stabilization (∆RE(σπ x ) and ∆RE(σπ y )) holds notable energetic significance.This observation underscores the existence of a synergistic CT (sCT) effect arising from the interplay between σ donation and π backdonation.The sCT effect does not straightforwardly align with the additive and separate consideration of σ donation and π backdonation, and it cannot be solely explained through ligand-field arguments.Notably, in the case of model III, which most closely mimics a real P450, ∆RE(σπ x ) and ∆RE(σπ y ) display the highest magnitude among all models, indicating that the proximal and equatorial ligands play a prominent role in promoting σ-π resonance.It should also be noted that there are two distinct sets of σ-π resonance interactions (σ-π x and σ-π y ) due to the presence of two π-backdonation orbital-interaction pairs (Figure 1).Comparatively, many other ligands, unlike CO, can only form a single π-backdonation orbital-interaction pair, resulting in just one set of σ-π resonance.This inherent disparity in sCT leads to a relatively heightened level of energetic stabilization in the ferrous heme-CO complex, distinguishing it from other complexes in terms of its bonding characteristics.In VBSCF calculations, the ∆RE(πxπy) values exhibit near-zero values across all models.The data indicate that the contribution of π-π resonance to the overall energetic stabilization is notably smaller than that of σ-π resonance.Nevertheless, in BOVB calculations, there is a significant enhancement of π-π resonance (∆RE(πxπy) = 4.03 kcal/mol in model III) compared to the VBSCF value (0.47 kcal/mol).This enhanced π-π resonance contributes to the additional stabilization in ferrous heme-CO bonding.This is especially relevant as many other ligands possess only one π-backdonation orbital-interaction pair and therefore cannot establish any π-π resonance.Consequently, the sCT effect resulting from π-π resonance also contributes to the additional energetic stabilization of ferrous In VBSCF calculations, the ∆RE(π x π y ) values exhibit near-zero values across all models.The data indicate that the contribution of π-π resonance to the overall energetic stabilization is notably smaller than that of σ-π resonance.Nevertheless, in BOVB calculations, there is a significant enhancement of π-π resonance (∆RE(π x π y ) = 4.03 kcal/mol in model III) to the VBSCF value (0.47 kcal/mol).This enhanced π-π resonance contributes to the additional stabilization in ferrous heme-CO bonding.This is especially relevant as many other ligands possess only one π-backdonation orbital-interaction pair and therefore cannot establish any π-π resonance.Consequently, the sCT effect resulting from π-π resonance also contributes to the additional energetic stabilization of ferrous heme-CO bonding.For comparison, we also computed the ∆RE(π x π y ) value for acetylene, yielding values of 13.41 kcal/mol (VBSCF) and 13.63 kcal/mol (BOVB) (Table S2).Thus, the energetic significance of π-π resonance is more pronounced in the case of acetylene.This finding also suggests that resonance between different bonding modes may be a widespread phenomenon extending beyond coordination complexes.
Scheme 3 provides a schematic representation of the sCT effect, a phenomenon uncovered in this study through ab initio VB calculations.sCT involves the contribution of bidirectional CT mechanisms.Firstly, resonance occurs between the VB structures within each of the σ and π frameworks, leading to pure CT stabilization.Secondly, the VB structures from distinct CT frameworks synergistically interact, conferring an additional stabilizing effect.The stabilizing σ-π resonance effect is particularly pronounced in ferrous heme-CO bonding, where two pairs of σ-π resonance can be established.The sCT effect is also feasible between two VB structures associated with π backdonation, although it was observed to possess less energetic significance compared to σ-π resonance.
Molecules 2024, 29, x FOR PEER REVIEW 9 of 13 Scheme 3. Schematic illustration of the sCT effect, arising from the resonance interaction between different sets of VB structures.

Materials and Methods
In our previous study [10], we optimized the geometry of a P450 ferrous heme-CO complex model containing porphine and HS -ligands at the B3LYP-D3BJ/def2-TZVP(6D,10F) level of theory [30][31][32][33][34][35][36].To better align with the conventional notation of dorbitals in chemistry such as  and  , we reorientated the optimized geometry such that Fe was placed at the origin, with Fe and the C atom of CO aligned along the z axis and one of the equatorial N atoms situated within the xz plane.For this purpose, we selected the N atom with the smallest magnitude of the N-Fe-S-H dihedral angle.The resultant model is referred to as model 0. However, applying VB calculations to model 0 was deemed computationally intensive.Therefore, we built three simplified models (I-III) based on model 0, as outlined in Scheme 1.These models retained the essential influence of proximal and equatorial ligands on the iron center while being computationally more manageable.The simplest model, model I, lacked any proximal or equatorial ligands, model II included a proximal HS -ligand, and model III featured both proximal and equatorial ligands (HS -and NH3).The atomic positions for models I and II were directly derived from model 0. We followed a similar procedure for constructing model III, but we optimized the positions of the hydrogen atoms within the NH3 ligands at the B3LYP-D3BJ/def2-TZVP(6D,10F) level using Gaussian 16 [37], while keeping all other atoms fixed.
Subsequently, we performed ab initio VB calculations for models I-III using the VBSCF method [38,39], employing the 6-31G*(6D,10F) basis set.VBSCF calculations allow for the optimization of both orbitals and VB structure coefficients {CK} within the total wave function Ψ: Scheme 3. Schematic illustration of the sCT effect, arising from the resonance interaction between different sets of VB structures.

Materials and Methods
In our previous study [10], we optimized the geometry of a P450 ferrous heme-CO complex model containing porphine and HS -ligands at the B3LYP-D3BJ/def2-TZVP(6D,10F) level of theory [30][31][32][33][34][35][36].To better align with the conventional notation of d-orbitals in chemistry such as d z 2 and d xz , we reorientated the optimized geometry such that Fe was placed at the origin, with Fe and the C atom of CO aligned along the z axis and one of the equatorial N atoms situated within the xz plane.For this purpose, we selected the N atom with the smallest magnitude of the N-Fe-S-H dihedral angle.The resultant model is referred to as model 0. However, applying VB calculations to model 0 was deemed computationally intensive.Therefore, we built three simplified models (I-III) based on model 0, as outlined in Scheme 1.These models retained the essential influence of proximal and equatorial ligands on the iron center while being computationally more manageable.The simplest model, model I, lacked any proximal or equatorial ligands, model II included a proximal HS -ligand, and model III featured both proximal and equatorial ligands (HS - and NH 3 ).The atomic positions for models I and II were directly derived from model 0. We followed a similar procedure for constructing model III, but we optimized the positions of the hydrogen atoms within the NH 3 ligands at the B3LYP-D3BJ/def2-TZVP(6D,10F) level using Gaussian 16 [37], while keeping all other atoms fixed.mixing of different structures for different CT routes.It should be noted that the contribution of the two-electron-transferred state (e.g., the third structure in Scheme 3 in the case of VB(σ)), was typically very small.
To effectively incorporate dynamical electron correlation into our wave functions, we extended our calculations with the BOVB method [46][47][48], building upon the VBSCF results.Again, we concisely described the non-VB electrons using a single set of MOs, focusing our explicit VB treatment on the active electrons.In this approach, VB orbitals are optimized separately for different VB structures.To enhance computational efficiency, our BOVB calculations selectively incorporated VB structures deemed relevant to either σ donation or π backdonation.In the VB(σ), VB(π x ), and VB(π y ) calculations with BOVB, we considered three distinct VB structures within their respective σ or π routes (see Scheme 3).Meanwhile, the VB(σπ x ), VB(σπ y ), and VB(π x π y ) calculations comprised 9 (=3 × 3) VB structures each, and VB(σπ x π y ) calculations encompassed a total of 27 (=3 × 3 × 3) VB structures.

Conclusions
In summary, our investigation into the nature of ferrous heme-CO bonding within cytochrome P450, employing ab initio VB calculations, has yielded valuable quantum mechanical insights.Firstly, our results underscore the dominance of π backdonation in the P450 ferrous heme-CO complex, primarily attributed to the presence of proximal and equatorial ligands.Furthermore, our VB findings illuminate the pivotal role played by σ-π resonance in stabilizing coordination bonding, which gives rise to the emergence of the stabilizing sCT effect.Notably, the presence of two pairs of σ-π resonance in ferrous heme-CO bonding, in contrast to the presence of just one in many other ligands, reinforces the bonding through augmented sCT.Additionally, one pair of π-π resonance also contributes to the bonding.By highlighting sCT as an adhesive force in coordination bonding, this study introduces a fresh perspective on bonding and opens new avenues of inquiry into coordination chemistry-a foundational field with widespread implications across various subdisciplines of chemistry.

Figure 1 .
Figure 1.Schematic illustration of one σ-donation and two π-backdonation orbital-interaction pairs for CT in ferrous heme-CO bonding.

Figure 1 .
Figure 1.Schematic illustration of one σ-donation and two π-backdonation orbital-interaction pairs for CT in ferrous heme-CO bonding.

Scheme 1 .
Scheme 1. Models I-III for VB analysis, constructed utilizing the previously optimized geometry of a P450 model[10].

Scheme 2 .
Scheme 2. Schematic drawings of representative VB structures.Black, red, and blue colors are used for electrons involved in σ donation, π backdonation in the xz plane, and π backdonation in the yz plane, respectively.

Scheme 2 .
Scheme 2. Schematic drawings of representative VB structures.Black, red, and blue colors are used for electrons involved in σ donation, π backdonation in the xz plane, and π backdonation in the yz plane, respectively.

Figure 3 .
Figure 3.The energy levels of the ten MOs around the frontier orbitals (B3LYP-D3BJ/def2-TZVP(6D,10F)) for models I-III without the CO ligand.The bars for occupied and unoccupied MOs are colored red and blue, respectively.

Figure 3 .
Figure 3.The energy levels of the ten MOs around the frontier orbitals (B3LYP-D3BJ/def2-TZVP(6D,10F)) for models I-III without the CO ligand.The bars for occupied and unoccupied MOs are colored red and blue, respectively.
a The total numbers of VB structures in VB(all) calculations with VBSCF and BOVB are 175 and 27, respectively.donationand the enhancement of π backdonation are even more pronounced in model III, where W2-W4 exhibit similar values (8.59, 8.52, and 8.53%, respectively), underscoring the role of equatorial ligands in promoting π backdonation.Another intriguing finding for model III is the relatively large values of W5 and W6 (9.27%), indicating that states concurrently involving one Heitler-London (HL) bonding pair for σ donation and one for π backdonation emerge as significant contributors to the overall wave function.Breathing orbital valence bond (BOVB) calculations produced results qualitatively similar to those obtained from VBSCF calculations (Table
a The total numbers of VB structures in VB(all) calculations with VBSCF and BOVB are 175 and 27, respectively.

Table 2 .
Effective weights (in %) of σ and π characters, determined via the VBSCF method.Values in parentheses are from BOVB calculations.

Table 3 .
RE values (in kcal/mol) for various states, obtained from (a) VBSCF and (b) BOVB calculations.